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Abstract. The present study deals with the problem of 
deriving the coefficients of the fundamental plane (FP) of 
early- type galaxies. 

We introduce statistical models of the FP and rela- 
tive fitting methods: the MIST (Measurement errors and 
Intrinsic Scatter Three dimensional) fits|^. The MIST fits 
account for the measurement errors on the variables and 
their correlations as well as for the intrinsic scatter. We 
show that the lack of a model of the intrinsic scatter of 
the FP is the origin of the systematic differences between 
the various fitting methods. 

We also address the problem of estimating the uncer- 
tainties of the FP coefficients and determine a simple re- 
lation between the sample size and the expected accuracy 
of the coefficients. 

The present study leads to define a 'minimum sample 
size' for a correct estimate of the uncertainties. For iV<30, 
both theoretical formulae and re-sampling techniques, like 
the bootstrap, do not give reliable estimates. 

The question of the 'universality' of the FP is ad- 
dressed by applying the MIST fits to ten samples of cluster 
galaxies. The FP slopes are actually consistent between 
the different samples, but, due to the large uncertainties, 
they could also hide significant systematic differences. 

The feasibility of the measurement of the possible vari- 
ations of the FP slopes as a function of redshift is also 
proved. 

Key words: Galaxies: fundamental parameters - Meth- 
ods: statistical - Galaxies: clusters: general. 



1. Introduction. 

The fundamental plane (FP) is a bivariate relation be- 
tween observed global properties of early-type (E) galax- 
ies like the effective radius r^ , the mean surface brightness 



/ie, and the central velocity dispersion ctq (e.g. Dressier et 
al. |1987|, Djorgovski & Davis |l987l) of the fornfl: 



log To = a log o-Q + bfic 



(1) 



The main features of the FP are its small scatter 10% 
in Te, of the same order of the measurement errors) and 
the 'tilt' of its slopes with respect to the prediction of the 
virial theorem for a family of homologous systems with 
constant M/L ratio. 

It is generally believed that the 'tilt' carries informa- 
tion on the nature of Es, although no convincing inter- 
pretation has been found that can be reconciled with the 
small scatter (see e.g. Ciotti et al. 1996, Rcnzini & Ciotti 
1993). It is now clear that part of the 'tilt' is due to the 
intrinsic non-homology, both structural and dynamical, 
of the E family (see e.g. Busarello et al. 1997, Graham 
& CoUess 1996, Prugniel & Simien [1997 and references 
therein). Stellar population effects (e.g. systematic differ- 
ences in age/metallicity and interstellar matter content) 
should account for the remaining tilt (see e.g. Mobasher et 
al. 199£). In that case, the tilt should decrease when mov- 



ing to near-infrared wavelengths, where these effects are 
minimized: some evidences exist that this is actually the 



case (e.g. Pahre et al. 199S, Scodeggio et al. 1998 and refer- 
ences therein) . The wavelength dependence of the FP en- 
ters also in the comparisons between FP determinations at 
different redshifts, where different rest- frame wavebands 
are sampled. 

The FP has been soon recognized as a precise tool for 
distance determinations (e.g. Dressier et al. 1987, Hudson 
et al. 1997) and, more recently, for constraining cosmo- 
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logical parameters and for studying galaxy evolution (van 
Dokkum et al. 1996, Bender et al. 1998| , J0rgensen et al. 
1999, Kelson et al. 200C, and references therein). 

All the above applications are based on comparisons 
between different determinations of the FP: e.g. different 
samples, different wavebands, and different redshifts. The 
method used to derive the FP coefRcients and the relative 



This representation of the FP is adopted throughout the 
paper. 
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uncertainties play therefore a central role for a proper use 
of the FP. 

The determination of the FP requires considerable ob- 
servational and data analysis efforts (see e.g. Ziegler & 
Bender 1997 and Ziegler et al. 1999), it is thus desirable 
that a similar effort be made for the derivation of the co- 
efficients of the FP and for an accurate estimate of the 
relative uncertainties. 

There is still no agreement in the scientific community 
about the fitting method to adopt. This is partly due to 
the numerous applications for which the FP is derived (see 
above), but it also originates from two important points: 

(1) the measurement errors on the FP variables are com- 
parable (see e.g. Smith et al. 1997, hereafter SLH97) and 
correlated (see J0rgensen et al. 1995a, hereafter JFK95a), 

(2) the scatter around the plane is not completely ac- 
counted for by the measurement uncertaint ies bu t has also 
an intrinsic origin (see e.g. J0rgensen et al. 1996 , hereafter 
JFK96). 

In the present work we address the problems relative 
to: (I) the choice of the procedure to derive the FP coeffi- 
cients (Sect. 2) and (II) techniques to estimate the corre- 
sponding uncertainties (Sect. 3). 

Starting from the bi-dimensional models introduced by 
Akritas & Bershady ( |1996| , hereafter AkB96), we propose 
in Sect. 2.1 statistical models of the FP that account for 
the various sources of scatter on the variables and derive 
the relative fitting procedures: the MIST {Measurement 
errors and Intrinsic Scatter Three dimensional) fits. 

We then discuss in Sect. 2.2 the capability of the fitting 
procedures to derive the coefficients of the FP. 

Section 3 deals with the uncertainties of the coefficients 
of the FP. The simulation algorithm used for the present 
analysis is described in Sect. 3.1. In Sect. 3.2 the different 
methods to estimate the uncertainties are analyzed and 
compared to the results of the simulations. 

The results and the techniques developed in Sect. 2, 
3, are applied in Sect. 4 to the comparison of the FPs of 
different clusters of galaxies. 

2. Deriving the FP. 

The choice of the best fitting procedure in the study of 
astronomical data is generally a not trivial question: the 
appropriate method should be always suggested by the 
scientific problem to be analyzed (see the discussion in 
Isobe et al. |1990[ , hereafter IFA90). 

The derivation of the FP is usually based on least 
squares methods: the ordinary least squares (OLS) fits 
(e.g. Lucey et al. 1991, Kjaergaard et al. 1993, Hudson et 
al. 1997), in which the root mean square (rms) of residuals 
relative to one of the variables is minimized, the orthogo- 



the arithmetic mean of the OLS coefficients (e.g. Faber 
et al. 1987). To reduce the effect of outliers, the fits are 



often performed by robust procedures, in which the sum 
of absolute residuals is minimized (see JFK96). 

Fitting methods based on multivariate analysis tech- 
niques, like the principal component analysis (e.g. Ben- 
der et al. 1992| ), have been sometimes adopted. By their 
nature, however, these methods are not suitable for the 
determination of best-fit coefficients, so that they will be 
not considered in the following. 

The most natural interpretation of the FP is in terms 
of a 'mean relation' between global quantities of E galax- 
ies with respect to which they scatter in the space of the 
observed properties. Drawing hint from AkB96, we will 
take this interpretation into mathematical terms, by in- 
troducing a statistical model for the FP (Sect. 2.1). The 
corresponding fitting procedures will be also derived. 

We will then discuss (Sect. 2.2) the problem of recov- 
ering the mean relation and illustrate the origin of the 
dependence of the FP coefficients on the fitting method. 

2.1. A statistical model: the MIST fits. 

Let us introduce three random variables {Xi}^^-^ ^ ^, de- 
scribing the distribution of global quantities of E galaxies. 
We will assume that {Xi} verify the following identity: 



^3 ^ aXi + fiX2 + 7 



(2) 



nal least squares (ORLS) fit (e.g. Busarello et al. 1997), in 



which it is minimized the rms of residuals perpendicular 
to the plane, and other methods derived from th e prev ious 
ones, as the bisector fit (e.g. Graham & CoUess 1996 ) and 



where a and (3 are the slopes and 7 is the zero point of 
the mean relation. 

In the following, we will use capitals to indicate ran- 
dom variables (RV) and the corresponding small letters for 
their outputs. Given two RVs A and _B, we will indicate 
with C{A, B) their covariance, with Y{A) and V(i3) their 
variances and with E(A) and E(i3) their expected values. 
The estimators of a given quantity (i.e. the RVs used to 
approximate the quantity we look for) will be marked with 
a caret (e.g. d, (3 and 7). An estimator will be said 'unbi- 
ased' if its expected value coincides with the quantity to 
be estimated. 

In modeling the scatter around the plane defined by 
Eq. (H) , we have not only to consider measurement errors 
on the variables but also an intrinsic dispersion. 

To this aim, let us introduce two sets of RVs, 
2 3 and {'I'mili^i 2 3' with zero expected values, 
that describe respectively the intrinsic dispersion and the 
scatter due to measurement errors, and let us consider the 
following relations: 

Y,^X, + + (3) 

where {^1}^^]^ 2 3 T(Sfs that describe the dis- 

tribution of the observed quantities of E galaxies (e.g. 
Yi =logro, Y2 = >3 = logCTo). 

The problem is then to determine unbiased estimators 
of a, (3 and 7. 
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From Eq. (||), one obtains the following identities: V(lj) = ^ . (jj^ _ Yj^ 



(15) 



E(X3) a • E{Xi) + 13 ■ E(X2) + 7 (4) 
C(X2,X3) = a ■C{XuX2)+ 13 ■C{X2,X2) = 

= a-C(Xi,X2)+/5-V(X2) (5) 
C(Xi,X3) = a-C(Xi,Xi)+/5-C(Xi,X2) = 

= a-V(Xi)+/5-C(Xi,X2) (6) 

where C(Xi,Xj) and V{Xk) are the components of the 
covariance matrix (CM) of {^i}- Assuming that {-'^"i}, 
{$mi} and {$si} are mutually not correlated (hereafter 
'hypothesis Hi') and making use of Eqs. (|) and (|), we 
can express the first and second order moments of {Xi\ as 
(see App. A for a straight demonstration and compare 
AkB96): 

nYi) = E(X,) (7) 
Y{Yi) = Y{Xi) + V($,n.) + V($s,) (8) 

C{Y,,Yj) = C(X„Xj)+C($,n.,$mj) + C($s„*s,) (9) 

Eliminating C{X,,Xj) and E(Xi) from Eqs. by 
means of Eqs. (0-^, we obtain a linear system of three 
equations in a, (3 and 7, whose solution is: 

a = • [C{Y2,Y^) - C(<i>^2, *m3) - C($s2, *S3)] ' 

• [C(yi, ^2) - C($,ni, *m2) - C($sl, $S2)] + 

- • [C(Yi,y3) - C($„u,$m3) - C(<fsl,*s3)] • 

■ [¥(^2) - v(a>,„2) - v(a>s2)] (10) 
/3 = A-i.[c(ri,y2)-c(<i>,ni,<fm2)-c($,i,<i>,2)]- 

• [C(yi, 1-3) - C(<f,ni, <i>,„3) - C($sl, $33)] + 

- • [C(r2, 1^3) - C(<i>„,2, <i'm3) - C($,2, *S3)] ' 

• [v(yi) -¥($„,! )-v(<i>si)] (11) 
7 = E(r3) - «E(yi) - /3E(r2) (12) 

where we put: 

A = [C(Yi,r2)-C($n,i,*m2)-C($sl,$s2)]' + 

- [V(Yi)-V($,nl)-V($,i)]- 

• [V(r2)-v(<i>,n2)-v(<i>.2)] (13) 

On increasing the sample size, the expressions obtained by 
substituting in the previous equations the terms G{Yi^ Yj) 
with their unbiased estimates C(yi,Yj) (hereafter d, (3 
and 7), will approximate more and more closely the 'true 
coefficients' a, (3 and 7. The quantities d, (3 and 7 will 
thus furnish unbiased, asymptotically normal, estimates 
of a, (3 and 7. In App. B the analytical formulae of the 
relative CM components are also derived. 
If {Yi} are normally distributed, the unbiased estimators 
for the moments of {Yi] are defined by the following usual 
expressions: 



G{Y„Y,) = ^ . ^ (y,, - Y,) ■ (yfe, - Y,) (16) 



i=l,N 



The fitting procedure based on Eqs. (|lOHl3|) and (yj-gj) 
will be indicated hereafter by the acronyms MIST, Mea- 
surement errors and Intrinsic Scatter Three dimensional 
fit. 



2.2. Deriving the coefficients of the FP. 



To discuss the application of Eqs. (10 -Ira) to the derivation 



(14) 



1=1, JV 



of the FP, it is important to remark the following points: 

TZi. The terms ¥($33) and V(<I>ni3) do not appear in 
Eqs. (p^-p^. The estimators of a, j3 and 7 are there- 
fore independent of the scatter along the dependent 
variable. 

TZ2. Setting equal to zero the other components of the 
CMs of J and in Eqs. (|lof-^ and using 

Eqs.(l^-Q), we obtain exactly the expressions of the 
OLSyg estimators, where the subscript ¥3 is used to 
indicate the dependent variable in the fit^. If, on the 
other hand, every variable is affected by dispersion, 
the quantities C(Yi,l2) would give a biased estimate 
of C(Xi, ^2) as shown by Eqs. (||, ||). To correct for 
the bias the CM components of {$mi} and {$si} are 
needed (see AkB96). 

Hz- As follows from TZi and 7?,2, to recover the mean 
relation by Eqs. (|l^-|l|) we have altogether ten de- 
grees of freedom, given by five elements of the CMs of 
{$nii} and {$si} respectively. 

TZi. Since Eqs. ^have the same dependency on the three 
variables, the best-fit coefficients obtained by means 
of Eqs. (f7^-|7^ are independent of the choice of the 
dependent variable. 

TZ^ . New unbiased estimators of the coefficients in Eq. (j2|) 
can be defined by taking some average of the slopes of 
the three planes determined by Eqs. (pO|~p^) using each 
of the {Yi} as dependent variable. We c an defi ne a bi- 
sector plane (see also Graham & Colless 1996| ), whose 
slopes are given by the vectorial sum of the normal 
vectors to the three planes obtained by Eqs. (|lo|-p^), 
and whose zero point is obtained by Eq. (^2|). The ad- 
vantage of the 'bisector' fit (hereafter MISTbis fit) 
consists in the larger effectiveness (i.e. smaller vari- 
ances for given sample size) of the relative estimators 
(see Sect. 3.2, and IFA90 for the two-dimensional case). 

TZq. There are often large systematic in- homogeneities be- 
tween different samples of data, due to differences in 
the procedures adopted to derive them (see e.g. Smith 

The OLS method can be thus regarded, in this respect, 
as a particular case of the MIST fit, in which the scatter is 
completely ascribed to the direction of the dependent variable. 
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et al. 1997). To account for such in-homogeneities, the 
least squares fits (ORLS and OLS) are often performed 
using robust estimators. The robust statistics can also 
be implemented to the MIST regression, by using ro- 
bust estimators for the moments of {Yi} in Eqs. (^0|- 

From the derivation of the equations in the previous 
Section, it follows that, in order to recover the coefficients 
of the mean relation, a, (3 and 7, some requests must be 
satisfied: 

7^1 : {Xi}, {$nij} and {$si} are not mutually correlated; 
7^2: estimates of the CM components of {$mi} are 
known; 

7^3: estimates of the CM components of {<i>si} are avail- 
able; 

7^4: unbiased estimates of the CM components of {Fi} are 
known. 

The hypothesis TLi is practically equivalent to assume 
that the measurement errors and the intrinsic scatter on 
the FP variables do not depend on the 'location' on the 
FP. This is not the case for the measurement uncertainties 
on the observed parameters and, a priori, it could be not 
true even for the intrinsic dispersion. However Eqs. 
[l3| ) continue to be valid by simply substituting the CM 
components of {$mi} and {$si} with their expected 
values as a function of the 'position' on the plane (see 
App. A and B for details). 

Concerning the measurement uncertainties (hypothe- 
sis 7^2), Eqs. (p^|-[l3|) can still be used by adopting aver- 
aged values of C($nij, $inj). To this aim, we suggest to 
proceed as follows. The quantities can be esti- 

mated as the square of the mean erroi^ on the parame- 
ters Ymi. For what concerns C(<I>nii, *I'mj), the only term 
that does not vanish, is C(logre, /Lte) due to the corre- 
lation between the uncertainties in log re and /^o- The 
quantities C(<I>ini, are never given in literature, so 

that one is forced to make some approximations. Since 
A log To ~ a-A/io with a ~ 0.3 (see JFK96), we can set 
C(logrc,//c) ~ a'V(/ic)- The uncertainties introduced by 
this approximation will be discussed in Sect. 4. 

Assumptions 7^3 and 7^4, cannot be satisfied because 
we do not have a physical model of the probability dis- 
tribution of E galaxies in the space of the observed quan- 
tities. It is thus necessary to introduce some simplifying 
assumptions: 

7i'3 : on the basis of TZi we can only assign all the intrinsic 

dispersion to the dependent variable. 
7i'4 : the only possible simple assumption is that {Yi} are 

normally distributed (see Eqs. |r^-|l6|). 

Both these assumptions introduce a 'bias' in the estimate 
of a, (3 and 7, so that the coefficients a, b and c in Eq. 

Notice that the mean errors are often the only significant 
available estimates. 



0-002 0-003 0-004 005 OOfi 



0-OOV 0-008 0-009 



Fig. 1. Coefficient a of FP simulations against the logCTo 
variance, fT?„„^ . The size of the simulations is TV = 3000. 

' log cro 

The symbols correspond to different fitting methods as 
shown in the lower left. 



do not necessarily coincide with a, (3 and 7 . Moreover, 
because of 7^'3, the estimates obtained by Eqs. (p^-p^ 
with a different choice of the dependent variable do not 
correspond to the same statistical model and so do not 
define the same plane. 

In order to better understand the bias introduced by 
7^'3 and 7^'4 , we used numerical simulations of the FP. 
The simulations were constructed giving a scatter around 
the plane in the logCTo direction by a normal RV (see 
Sect. 3.1 for details). The variance of logCTo has been var- 
ied in the range of values obtained for the log ctq scatter of 
different samples of galaxies (see Sect. 4.1). The FP coeffi- 
cients have been derived for each simulation by the MIST 
fits, with assumptions 7^'3 and 7^'4, using each of {Yi} as 
the dependent variable (hereafter MIST-j^^ fits) , and by the 
MISTbis- We used the estimates of C($mi, $mj) obtained 
by typical values of the uncertainties on FP parameters 
(0.03 in logTe, 0.1 in and 0.03 in logCTp). 

In Fig. |l|, we plot the coefficient a against CTj^^ (sim- 
ilar results are obtained for h and c) . The values of the FP 
coefficients obtained by the various fitting methods turn 
out to be systematically different, by an amount that in- 
creases with the scatter. An estimate of this difference for 
the MIST fits will be given in Sect. 4 by comparing the 
FPs of different clusters of galaxies. 

The assumption TL' 4^ introduces likely the same 
amount of bias in the various fits. To illustrate this point, 
let us consider an example in which the CM components 
of the intrinsic dispersion are 'known'. We will assume 

V($„Uog.o)+V(*«logao) = V($„Uogre)+V($slogrJ = 

0.0035 (a scatter of ^ 14% in re and cto)^ and that the 

^ This amounts to a 15% scatter in the ctq direction, i.e. to 
the mean value of the FP scatter of clusters, see Sect. 4. 
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Table 1. Values of a, b and c of a FP simulation with intrinsic dispersion 'known' (see text). The sample size is 
N — 3000. Column 1: FP coefficients a, b and c. Columns 2, 3, 4 and 5: values obtained by the different MIST fits 
(seesubscripts in the top row) with corresponding uncertainties {la intervals). 





logic 




logcro 


BIS 


a 


1.26 ±0.02 


1.28 ±0.02 


1.28 ±0.01 


1.28 ±0.01 


b 


0.310 ±0.002 


0.308 ± 0.002 


0.309 ± 0.002 


0.309 ± 0.002 


c 


-8.27 ±0.07 


-8.23 ±0.07 


-8.20 ±0.07 


-8.24 ±0.06 



other CM terms of {$mi} and {$si} vanish. We now 
create a FP simulation, and derive a, b and c by the MIST 
fits, using the assigned values of the CM components of 
{$si} and {$mi}, and making therefore onZy the hypoth- 
esis T-l'i- The results of the fits are shown in Tab. |l|: the 
FP coefficients are practically independent of the fitting 
method. 

We then conclude that the differences of Fig. |l| are due 
only to the assumption Ti.'^: the lack of a model for the 
intrinsic dispersion of the FP variables is thus at the ori- 
gin of the dependency of the FP coefficients on the fitting 
procedure. The larger is the scatter around the plane, the 
larger will be the bias due to the fitting method. 



3. Uncertainties on the FP coefficients. 

Two kinds of methods are usually adopted to estimate 
the uncertainties on fit coefficients: theoretical methods 
and re-sampling techniques. 

The 'theoretical uncertainties' are obtained from the 
analytical expression of the variances of the estimators 
(e.g. IFA90 and Feigelson & Babu 1992). By their nature, 



these estimates are valid only asymptotically, i.e. for large 
sample sizes. 

When analytical formulae are not available, or when 
the sample is small, re-sampling procedures are adopted. 
The statistics of interest is calculated for various 'pseudo- 
samples' drawn from the original data set. The uncer- 
tainties are then derived from the distribution of pseudo- 
values. The main re-sampling procedures are known as 
'jackknife' (see Quenouille 194£ and Tukey 1958|) an d 
'bootstrap' (see Efron |l979| and Efron & Tibshirani 19861 ). 
In the jackknife, one point is extracted in sequence from 
the original data set, so that a number of pseudo-samples 
equal to the sample size is constructed. In the bootstrap, 
random samples are drawn by replacement from the actual 
data set. 

Figure ^ shows the relative uncertainties on the FP 
coefficients a and &, as drawn from literature, as a function 
of the sample size N. The large scatter of the values of the 
uncertainties is indicative of the inconsistency between the 
methods used in the different works to estimate the errors. 




Fig. 2. Relative uncertainties Sa/a and 6b/b on FP slopes 
against the sample size (see text). The values are drawn 
from literature. Sources of the uncertainties and relative 
s ampl e sizes (within square brackets): Busarello et al. 
( |1997| ) [40]; de Carvalho & Djorgovski ( |l992[ ) [31 4 55 
58]; D'On ofrio et al. ( |l997D [25]; Dressier et al. ( |1987| ) [43]; 
Graham ( |1998D [7 18 25 25]; Graham & CoUes s (|l996|) [26]; 
Hudson et al. ( |1997D [352]; J0rgensen et al. ( |1996D 79 22 
24 14 10 14 29 8 19 226 109 109 109 1 09 41 41]; Mobasher 
et al. ( |1999| ) [48]; Pahre et al. ( |1995[ ) [8 10 12 13 16 59]; 
Pahre et al. ( |1998D [251 60 16 23 15 14 15 17 11 19 11 27] ; 
Prugniel & Simien ( il994| ) [102]; Prug niel & Simien ( |1996| ) 
[231 39 21]; Recill as-Cruz et al. ( |1990D [29 29 29] ; Rcci Uas- 
Cruz et al. ( |l99lD [22 22 22]; Scodeggio et al. (|1998D [38 
41 54 75 73 29]. 



In the following we will analyze the performance of 
the methods to estimate the FP uncertainties with special 
regard to the role of the sample size. The analysis will be 
performed by numerical simulations, which are described 
in the next section. 
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Fig. 3. Comparison of one of the FP simulations with the template sample, (a) Simulation and the Coma photometric 
sample in the plane (/io,logrc). (b) Distribution in the plane {fic, logao). (c) Distribution in the plane (logrclogcro)- 
(d) Face-on view of the FP. (e), (f) Long and short edge-on views of the FP. The plotted simulation has the same size 
as the Coma sample. 



3.1. The simulation algorithm. 

The simulations consist in distributions of points ex- 
tracted from a common 'parent distribution'. To derive 
the parent population we based on the distribution in 
the parametric space of the galaxies in the Coma cluster 
(hereafter the 'template' sample). This choice is manda- 
tory since Coma is the only cluster with FP parameters 
available for a large number of galaxies. 

The simulation algorithm consists in the following. 

At first, we derived a parent distribution of E galaxies 
in the plane (/ic,log To), using the photometric data in 
the Gunn r band by JFK95a for 146 galaxies in Coma. 
The sample is complete out to Gunn r ~ 15.5 mag. 
The distribution with respect to /Xe was described by 
a normal RV. The fic interval of the template was 
binned, and the mean value (MV) and standard devi- 
ation (SD) of the logTc distribution derived for each 
bin. The MVs and the SDs were then fitted with re- 



spect to the central values of /ic by polynomials of 
suitable order, and the best-fit curves were used to in- 
terpolate a value (and a scatter) of log to each value 
of Re- 
using (To by J0rgensen et al. ( |l9951j , hereafter JFK95b) 
for 75 galaxies of the photometric sample (see Sect. 4), 
we determined the FP coefficients and the root mean 
square (rms) of logcro residuals, siogao, by the 
OLSiog^o fit. We found a ~ 1.50, b ~ 0.328, c ~ -9.1 
and siogaa — 0.048. These quantities were used to as- 
sign to the points a value and a scatter with respect 
to the variable logcro. 



One of the simulated samples is compared to the tem- 
plate in Fig. 1^, to show how the simulations resemble very 
well the distribution of Coma galaxies. 
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3.2. Estimating the uncertainties. 

We derived the 'true uncertainties' on the FP coefficients 
as a function of the sample size in the following way. 

FP simulations of fixed size A'' were constructed and 
the coefRcients a, b and c determined for each sample. 
The 'true uncertainties', Sa, Sb and 6c were determined 
from the distributions of a, b and c by using la standard 
intervals. To obtain estimates independent of the 'fitting 
scale', we divided Sa, 6b and 6c by the 'true' coefficients 
of the FP, at , bt and Ct , derived from a simulation of large 
size. In Fig.^ we plot 6a/at, 6b/bt and 6c/ct against the 
logarithm of the sample size N for the various MIST fits. 
To allow a direct comparison with Fig. ^ the same range 
of N has been plotted. 

It is apparent that the true intervals depend on the fit- 
ting method and that the most effective fit (i.e. lower val- 
ues of 6a/a and 6b/b for fixed N) is the MISTbis method, 
in agreement with what found by IFA90 for the bisector 
line. The largest variances are obtained for the MISTiogo-o 
fit. It is also worth to notice that, by changing the scat- 
ter along the ctq direction in the simulations, the curves 
shown in the figure undergo just a translation in the y di- 
rection, without any change in shape. In the case shown in 
Fig. H we adopted the MV of the rms of log ctq residuals 
of various clusters of galaxies (see Sect. 4). 

Figure ^ can be used as a ready tool to state the num- 
ber of galaxies necessary to achieve a given accuracy on 
the FP. 

Concerning the zero point of the FP, it is iiriportant 
to remark that the uncertainties plotted in Fig.Q do not 
represent the estimates of usual interest. For all the appli- 
cations of the FP zero point (i.e. distance determinations, 
constraining of cosmological and evolution parameters), 
the uncertainties on c are derived with the hypothesis that 
the FP slopes are exactly determined. The errors on c are 
then given as 1/ ^/N multiplied by the rms of the de- 
pendent variable residuals. Such estimates are generally 
smaller than 6c/ Ct. 

For this reason, in the following we will focus our anal- 
ysis on the uncertainties of the FP slopes. 

Although the comparison of Figs. |^and ^ does not 
show an evident disagreement, for logA^<1.8 {N<60) the 
values reported in literature appear almost as a scatter 
diagram. Starting from this remark, we now analyze the 
performances of the different methods used to estimate 
the uncertainties. 



3.2.1. Theoretical methods. 

Although statistics allows to prove the asymptotic validity 
of variance estimators, it does not furnish an estimate of 
the 'minimum sample size' for the theoretical formulae to 
be valid. Since such an estimate will generally depend on 
the 'shape' of the parent population, it should be obtained 
each time by using simulation methods. 




Fig. 4. 'True relative uncertainties' on the FP coefficients 
against the logarithm of the sample size N. The various 
MIST fits are plotted by different symbols as shown in the 
upper right-hand corner. The MISTiog^, and MIST^^ 
fits are completely consistent, so that for clarity only the 
first one is plotted. The axes scales are the same as in 
Fig. |. 

To test the performance of theoretical variance estima- 
tors for the FP coefficients we apply the results of Sect. 2. 
FP simulations of fixed size N were constructed and the- 
oretical relative uncertainties {la standard intervals) on 
a and b derived for each sample by Eqs. (B17, |l^ - p^ . 
In Fig. H, the MVs of the theoretical estimates, 6a/ a and 
(5&/&, with corresponding 'error bars', are plotted against 
the logarithm of N. The error bars were obtained by con- 
necting the 5th and 95th percentiles^ of the 6a/a and 6b/b 
distributions. For comparison, we also plot the 'true un- 
certainties' on the slopes as derived from the simulations. 
As expected, the larger the sample size, the better it is 
the agreement between theoretical and actual values. For 
sample sizes smaller than N^i^ « 30 (logA^min = 1.5), we 
see however that theoretical formulae (I) largely scatter 
(up to 40% for a and 15% for b) and (II) increasingly un- 
deresti mate th e 'true' values (see also IFA90 and Feigelson 
& Babu |1992[ ). 

To discuss in more detail the point II, we calculated for 
each sample the coefficient a, the theoretical uncertainty 
(5a, and the 'discrepancy' da = \a — at\/6a, where the ac- 
tual value, at, of the coefficient a was derived by fitting a 
simulation of large size. We then determined the fraction 
of simulations, /a, with da greater than a fixed value A. 
The calculation was iterated by varying the sample size. 

In Fig. H, we plot /a against logA^ for two different 
values of A. We chose A = 1.6 that, for a normal distribu- 
tion, defines a confidence level (CL) of 10%, and A = 2.6, 

® The pth percentile of a given distribution is the point that 
exceeds (lOO-p)ascending order (see e.g. Beers et al. 1990). 
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logN logN 



Fig. 5. Comparison of theoretic and true relative uncer- 
tainties on the FP slopes (see text). The different symbols 
are explained in the upper right-hand corner. We adopted 
the MISTiogr, fit. 



corresponding to a f % CL. If Eqs. (|T3 - H) worked well, 
on the average, for every value of N, the fraction of not- 
consistent samples would be independent of N and deter- 
mined by the CL corresponding to the value of A. Figure ^ 
shows that for N > N^nin ~ 30 (logA^niin ~ 1-5) /a is in 
good agreement (within ~5%) with the expected CLs. For 
sample sizes smaller than N^i^i the theoretical formulae 
do not furnish reliable estimates of the desired CLs: the 
curves of Fig. |^ increase steeply. The same result was ob- 
tained for the coefficient b, and by varying the MIST fit. 
We also found that the estimate of Ny^in is largely inde- 
pendent of the simulation parameters: the same estimate 
was obtained using the FP coefficients and the log uo dis- 
persion of the samples studied in Sect. 4. 

Figure ^ suggests that in order to obtain a given con- 
fidence level, the uncertainties on the FP coefficients of 
small samples should be estimated using a a interval de- 
pendent on TV. For instance, to obtain a 10% CL with 
TV = 12 (log = 1.1), an effective interval of 2.6a should 
be used. However, it turns out that the smaller the sample 
size with respect to N,^nim the stronger is the dependence 
of /a on the adopted fitting method and, what is more 
critical, on the simulation parameters. For instance, the 
fraction of not-consistent samples with = 10 varies by 
~ 6% by varying the logcro scatter in the simulations. 

We conclude that for A^<30 the theoretical formulae 
are not reliable. Although the desired confidence intervals 
can be roughly obtained by using effective, suitably tested, 
standard intervals, the individual estimates can be signif- 
icantly different, up to ~ 20% for a and ^ 15% for b, with 
respect to the true values. 



Fig. 6. Fraction of non-consistent simulations, /a, 
against the sample size N (see text). Circles and rect- 
angles correspond to different values of A as shown in the 
upper right-hand corner. The solid lines were obtained by 
interpolation. The dashed lines show the asymptotic val- 
ues of the curves. 



3.2.2. Re-sampling procedures. 

The hypothesis underlying the use of re-sampling methods 
is that the available data set furnishes a good approxima- 
tion to the parent population. The statistics of interest 
is calculated for various pseudo-samples drawn from the 
actual data set. If this 'sampling hypothesis' holds, the 
distribution of pseudo-values coincides with the 'true' one 
and the confidence intervals can be accurately estimated 
at the cost of some computing time (see e.g. Efron & Tib- 
shirani 



1986). However, the smaller the sample size, the 



larger is the probability that the actual sample gives a 
poor representation of the parent population. In order to 
derive a 'minimum sample size' for the re-sampling meth- 
ods to be reliable, numerical simulations have to be em- 
ployed. 

On the basis of the analysis of the previous Section, 
we studied the performance of re-sampling uncertainties 
on the FP coefficients by testing the use of the bootstrap 
method. The results that follow were found to be largely 
independent of the simulation parameters, of the actual 
FP coefficients, and of the fitting method. 

FP simulations of fixed size N were constructed. For 
each sample, we determined the FP slopes and the boot- 
strap uncertainties Sa and Sb by using la standard inter- 
vals with iVr = 2000 pseudo-samples (by further increas- 
ing TVr the following results do not change). In Fig. ^ 
the MVs of the bootstrap relative uncertainties, 6a /a and 
Sb/b, with corresponding error bars (see Sect. 3.2.1), are 
plotted against the logarithm of the sample size and com- 
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Fig. 7. Comparison of 'bootstrap' and 'true' relative un- 
certainties on the FP slopes (see text). Different symbols 
have been used as indicated in the upper right-hand cor- 
ner. We adopted the MISTiogr-, fit. The axes scales are 
the same of Fig. ^. 



Fig. 8. Fraction of non-consistent simulations, fp, against 
logA^ (see text). The plot is referred to the coefficient 
b. Circles and rectangles correspond to 100p% percentile 
intervals with p = 0.10 and p = 0.01 (see upper right-hand 
corner). The dashed lines show the desired CLs. The axes 
scales are the same as in Fig. |[ 



pared to the true values. The error bars are given by 10% 
percentile intervalsQ 

As shown in the figure, the bootstrap method gives 
a good measure of the average uncertainties, but a very 
poor, largely scattered, representation of the actual errors. 
Only for iV<15 (logiV<l.l - 1.2) the MVs of the re- 
sampling standard errors appear to differ from the true 
uncertainties. As a matter of fact, it turns out that this 
difference is due to the use of la standard intervals to es- 
timate the bootstrap confidence intervals. For small sam- 
ples, in fact, the distribution of pseudo- values is signif- 
icantly different from a normal one, so that the desired 
CLs must be obtained by non-parametric estimates^ (see 

To illustrate 



1987) 



Efron & Tibshirani 1986 and Efron 
this point, we tested the use of bootstrap percentile in- 
tervals proceeding as in the previous Section (Fig. ||). For 
each sample size, we derived the fraction of simulations, 
fp, that are not consistent with the 'true' FP slopes. The 
calculation was iterated by adopting different percentile 
intervals of the pseudo- values. To have a direct compari- 
son with Fig. H we chose 100 • p% intervals with p = 0.10 
and p = 0.01 respectively. In Fig. ^, the fraction of non- 
consistent samples is plotted against the logarithm of the 
sample size. It turns out that, on the average, the boot- 
strap allows to give very accurate estimates of the true 



A p% percentile interval is given by the difference between 
the p/2th and the (100 — p/2)th percentiles (see note |^ of the 
distribution of interest. 

® I.e. by confidence intervals independent of any assump- 
tion on the shape of the probability density of the statistics of 
interest. 



intervals. For every sample size, differences of only some 
percents are found with respect to the desired CLs. 

On the other hand, by looking at Fig. ^ we notice 
that for small samples the bootstrap uncertainties have 
a very large dispersion with respect to true values. The 
scatter varies from ^ 5% for N = Nj^in ~ 30 up to 60% 
for N « 10. Below N « Nj^in the FP parent population is 
poorly represented, so that the single bootstrap estimates 
can be very unsatisfactory. 

To have a comparison with the theoretical methods, 
we compared the width of the error bars shown in Figs. |^ 
and 0. In Fig. || we plot the difference of the error bars 
against the logarithm of the sample size. 

While for large samples theoretic and bootstrap un- 
certainties have a similar scatter, for N<Njnin ~ 30 
(log A^min<1.5) the bootstrap errors become increasingly 
less accurate. For N « 10, the scatter increases up to 
- 30% for a and to - 20% for b. 

The conclusions are thus the following. 
For large samples, both theoretic and bootstrap meth- 
ods give accurate estimates of the uncertainties. 

For A^<30, both methods give estimates that can differ 
significantly from the true values. The bootstrap is accu- 
rate on the average but the uncertainties have a very large 
scatter. The theoretical methods give values that are more 
precise, but systematically underestimated. 
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Fig. 9. Differences, ea{bo — th) and £b{bo — th), of the 
'scatters' of bootstrap and theoretic uncertainties on a 
and b, against logiV. 



Fig. 10. FP slopes, a and b, of the samples listed in 
Tab. ||. The ellipses define a 10% CL. The coefficients were 
derived by the MISTiog r„ fit. 



4. An application: the FP of cluster galaxies. 

The use of the FP as distance indicator is based on the 
hypothesis that its slopes and thickness are 'universal', 
i.e., independent of the sample of galaxies. 

The problem of the universality has been studied by 
JFK96 by deriving the FPs of ten clusters of galaxies. Al- 
though the authors find that the FP coefficients are no 
significantly different from cluster to cluster, they admit 
that the samples are too small to get an accurate compar- 
ison. 

Up to now, the FP has been derived with a significant 
number of galaxies only for Coma cluster {N = 81). Most 
of th e othe r cluster samples have -/V<2 — 30 (e.g. D ressie r 

119981 ). 



et al. 1987, JFK96, Hudson et al. 1997, Pahre et al 



On the basis of the discussion in Sects. 2 and 3, we try 
now to address the question of the universality by com- 
paring the FPs of different clusters. 

We chose the samples so that a) their size was as large 
as possible; b) the FP parameters of the different samples 
were homogeneous; c) each sample had photometric pa- 
rameters in the same waveband and d) the cluster were at 
2^0. The four prescriptions are needed to a), b) make 
significant comparisons, c) avoid waveband dependencies, 
and d) avoid galaxy evolution effects. The ten samples 
surviving the above criteria are listed in Tab. ^. 

The FP parameters, logcrp, log re and /io were drawn 
from literature. The photometry is in Gunn r for JFK data 
and in Kron-Cousins R band for the other samples. For 
the present study, the differences of the two wavebands 
can be completely neglected (see Pahre et al. 1998 and 
Scodeggio et al. 1998). Only galaxies wi th gg > 95kms^^ 
were considered (see e.g. Dressier et al. 1987). 



The FP coefficients were determined by the MISTiog r„ , 
MIST^^ , MISTiog and MISTbis fits (see Sect. 2.1), by 
taking into account the measurement errors on the vari- 
ables as described in Sect. 2.2. 

To compare the FP slopes, it should be taken into ac- 
count that the uncertainties on a and b can be correlated. 
By Eqs. B17, we obtained the theoretical estimates of the 
CM components of the slopes. These estimates were then 
used to derive the ellipse that defines a 10% CL for a nor- 
mal bi-variate. This is accounted for in Fig. |l^. 

The Fig. |l^ has the disadvantage to be not easily read- 
able. To obtain a more immediate description we com- 
pared separately the a and b values as derived from dif- 
ferent fits, as plotted in Fig. ^ The FP slopes of the 
MISTiog r„ and MISTiog ^0 fits are also listed in Tab. |[ 

All the samples (except one) have small size, N < 



30, so that, as discussed in Sect. 3.2, an accurate 



estimate of the uncertainties on a and 6 is not achievable. 
Approximate error bars were derived by theoretic formu- 
lae (Eqs. B17), with the prescriptions of Sect. 3.2.1 to 
obtain a 10% CL. Since half of the clusters have iV<20, 
effective theoretic intervals should in fact be more suitable 
than re-sampling estimates (see Sect. 3.2.2). 

Although the coefficients a and b are consistent for each 
pair of samples in every fitting procedure, we see that the 
error bars are very large, typically 40% — 50% of the a 
and b values. Moreover, they could not represent reliable 
estimates. 

It is also worth to notice that in the MISTiog , 
MIST^^ and MISTbis fits, the coefficient a of the largest 
sample (the Coma cluster) is systematically higher with 
respect to the other determinations. A weak correlation 
between a and the sample size seems also to be present. 
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Fig. 11. FP slopes of the clusters listed in Tab. g. The coefficients a and b of the different MIST fits are plotted on 
horizontal panels. The fitting methods are indicated shortly in the lower right-hand corners of the panels on the left. 
The sample size increases from the left to the right of each panel as shown on the x-axes of the lower panels. The 
error bars define a 10% CL. The values of a and b after the magnitude-completeness corrections are also plotted by 
asterisks. The dashed lines correspond to the weighted means of the slopes. 
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Table 2. Column 1 : cluster identification. Column 2 : CMB-frame redshift. Column 3 : limiting apparent magnitudes. 
The values were referred to the distance of Coma by using the redshifts of column 2. Column 4 : number of galaxies 
(E+SO). Column 5 : references from which the FP parameters are drawn. Column 6, 7 : FP slopes, a and 6, obtained 
by the MISTiogr-„ fit, with Icr standard intervals. Column 8, 9 : FP slopes and relative uncertainites {la standard 
intervals) as derived by the MISTiogo-o method. Column 10 : intrinsic dispersion of the FP, Siogo-o, projected along the 
log (To direction. The values were obtained by the MISTiog^o fit. 



Id. 


cz 

(kms^^) 


mi 


iV 


Ref. 


a 


6 


a 


6 


Slog o-Q 


S639 


6545 


14.9 


12 


JFK95a,b 


1.14 ±0.20 


0.306 ± 0.033 


1.49 ±0.25 


0.300 ± 0.043 


0.058 


J8 


9800 


13.8 


13 


SLH97 


0.77 ±0.24 


0.317 ±0.022 


1.38 ±0.26 


0.330 ± 0.040 


0.070 


A3381 


11472 


14.8 


14 


JFK95a,b 


1.10 ±0.13 


0.249 ± 0.028 


1.50 ±0.26 


0.232 ± 0.026 


0.057 


S753 


4421 


14.7 


16 


JFK95a,b 


0.89 ±0.16 


0.341 ±0.020 


1.21 ±0.16 


0.353 ±0.021 


0.061 


HYDRA 


4050 


15.1 


18 


JFK95a,b 


1.19 ±0.14 


0.320 ±0.028 


1.90 ±0.45 


0.316 ±0.033 


0.062 


A194 


5038 


14.7 


21 


JFK95a,b 


0.94 ±0.12 


0.348 ± 0.025 


1.26 ±0.13 


0.344 ± 0.025 


0.058 


PISCES 


4700 


15.5 


24 


SLH97 


1.16 ±0.09 


0.348 ±0.037 


1.52 ±0.31 


0.348 ± 0.040 


0.056 


A539 


8734 


14.9 


25 


JFK95a,b 


1.15 ±0.12 


0.317 ±0.014 


1.51 ±0.23 


0.327 ±0.018 


0.029 


PERSEUS 


4800 


15.1 


30 


SLH97 


1.32 ±0.18 


0.348 ± 0.014 


1.68 ±0.17 


0.370 ±0.017 


0.010 


COMA 


7200 


15.4 


75 


JFK95a,b 


1.34 ±0.07 


0.326 ±0.012 


1.53 ±0.08 


0.328 ±0.013 


0.015 



To address this point, we tried to correct the FP slopes 
for the different magnitude-completeness of the samples 



(see Giovanelli et al. 1997 and Scodeggio et al. 1998 for 
a wide discussion). At first, we constructed the complete- 
ness histograms of each sample. To this aim, the mag- 
nitude range was binned and the fraction of galaxies of 
each bin normalised to the corresponding fraction of the 
Coma photometric sample of JFK95a. This data set is 
complete in fact out to a magnitude higher than the mag- 
nitude limits of the other samples (see Tab. || and refer- 
ences therein). It turns out that the same results are also 
obtained by normalizing the fractions of galaxies through 
a gaussian model of the luminosity function (see Scodeg- 
gio et al. 1998). The histograms were then modeled by 



Fermi-Dirac distributions and incomplete FP simulations 
constructed through the modeled distributions. By fitting 
the simulations, we estimated the corrections on a and b. 

As shown in Fig. the corrections shift upwards the 
coefficient a and reduce the correlation with the sample 
size and the systematic difference of the Coma sample. It 
turns out, in fact, that these effects are a consequence of 
the different completeness of the samples with respect to 
the photometric parameters. This is shown by the results 
of the 'inverse' fit, with logcro as the dependent variable, 
that should be less sen sitive to the photometric complete- 
ness (see Hudson et al. 1997 , hereafter HLS97). By looking 
we see that in the MISTiogo-o 



at Fig. 11 



method there 

is no systematic difference between Coma and the other 
samples. 

It is also evident that a systematic difference is present 
between the FP coefficients derived from the various MIST 
fits. As discussed in Sect. 2.2, the 'fitting bias' is due to 



the lack of a model for the intrinsic scatter of the FP. 
For the same reason, only the projection of the intrinsic 
dispersion along some assigned direction, with respect to 
a given fitted plane, can be measured. 

For instance, we calculate the scatter projection on the 
logCTo variable for MISTiogo-o fit. For each sample we de- 
rived the rms of the logtio residuals. The amount of scat- 
ter due to the measurement errors was then subtracted 
in quadrature. The measurement error scatter projected 
on (To was found to be ~ 10%. By constructing FP sim- 
ulations with (To ^ 100, the 'intrinsic rms values' were 
corrected for the bias due to the log (To cut on the cluster 
samples. The correction was found to be very small (~ 2 
- 3%) and largely independent of the simulation parame- 
ters. 

The log ctq dispersions are shown in Tab. |. The mean 
values of the logcro intrinsic scatter amounts to ^ 0.048, 
that corresponds to ~ 10% in (Jo- 
in Tab. |, we show the weighted means of the coeffi- 
cients a and h for the various fitting methods. The means 
were calculated after the magnitude-completeness and the 
log (To corrections. The logcro bias on the slopes of the 
MISTiogcro fit was estimated from the simulations (see 
above). It amounts to ^ 5% for a and is completely neg- 
ligible for h. 

The difference between the various MIST fits amounts 
to ~ 15% for the coefficient a and to ~ 6% for h. The 
coefficient a varies from 1.25 to 1.45, and h from 0.33 to 
0.35. 

The bias that could be introduced by neglecting the 
measurement errors was found to be negligible for b 
(<1%). Concerning the coefficient a, the bias depends 
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Table 3. Weighted means of the 'corrected' FP slopes 
(see text). Column 1: fitting method. Column 2 and 3: 
mean values of a and b with corresponding uncertainties 
(Ict intervals). 





a 


b 


MISTiog,^ 


1.25 ±0.05 


0.329 ± 0.008 


MIST^„ 


1.27 ±0.05 


0.350 ± 0.007 


MISTiog^o 


1.45 ±0.06 


0.332 ± 0.009 


MISTbis 


1.32 ±0.04 


0.339 ± 0.007 



Coma 



Fig. 12. Comparison of the FP slopes obtained by the 
MIST fits with the values of JFK96 (see text). Symbols 
for the different MIST fits are as in Fig. 0. 



on the fitting procedure: it amounts to ^ 5% for the 



MIST 



JFK96 appear to be consistent with our estimates for each 
cluster. However, some systematic difference exists. The 
values of a obtained by the MISTiog and MIST^^ fits 
are systematically lower than those of JFK96. A small sys- 
tematic difference is also found for the MISTbis method. 
For the logCTo fit the MIST estimates are systematically 
higher. 

The values of b agree with those by JFK96 for the 
MISTiog , MIST^„ and the MISTbis methods. On 
the other hand, the values of the MIST^;^ fit are system- 
atically higher. 

The values of the FP slopes obtained by JFK96 for the 
whole cluster sample, a = 1.24±0.07 and b = 0.328±0.008, 
are generally consistent with the mean values reported in 
Tab.|[ However, systematic differences are particularly ev- 



log To 



and MISTp^ methods (in agreement with the ident for the value of a in the MISTiog ag fit and for the 



estimate given by JFK96), and to ~ 1% and ~ 3 - 4% 
for the MISTiog ^0 and MISTbis fits. Because of the 
correlation between the uncertainties on logro and ^o, it 
turns out that the MISTiog o-q method is quite insensitive 
to the uncertainties on the FP parameters. 

For what concerns the covariance term of the photo- 
metric parameters, as discussed in Sect. 2.2, it was es- 
timated by the relation C(logro,//c) ~ ct ■ V(/ie) with 
a ~ 0.30 (see JFK95a). By analyzing the correlation of 
the errors on log and fj,o for various photometric data set 
drawn from literature, we obtained values of a included 
between 0.26 and 0.34. Varying a in this range for the 
calculation of C(logro, /io), very small variations (< 1%) 
were obtained in the FP coefficients. The higher system- 
atic difference, 2 - 3%, was found for the coefficient a 
in the logcrp fit. 

The systematic dependence of the FP coefficients from 
the fitting method can also be seen by comparing the val- 
ues obtained by the MIST fits with the results of JFK96 
and HLS97. 

In Fig. [l^ the MIST estimates of a and b are com- 
pared to the values of JFK96, that adopted an ORLS (ro- 
bust) method. To this aim, we plot the differences of the 
FP slopes (JFK96 — ours) divided by the relative Icr un- 
certaintities added in quadrature. The values obtained by 



value of b in the MIST^^ method. 

HLS97 derive the FP slopes for a sample of seven clus- 
ters of galaxies. By adopting the OLSiogo-Q method on all 
the clusters simultaneously, they obtain a — 1.383 ± 0.04 
and b = 0.326 ± 0.011. These values are consistent with 
the mean values of the MISTiog^,, fit (see Tab. 1). 
The slopes of the individual clusters of HLS97 are derived 
by a bi-dimensional fit, adopting as independent variable 
the combination of log and fi^ obtained from the global 
fit. A direct comparison of our estimates with the individ- 
ual results of HLS97 is thus not possible. 

Since the ORLS and OLS methods do not account for 
the measurement errors, and due to differences in com- 
pleteness and selection, a detailed explanation of the ori- 
gin of the above discrepancies is not achievable. 



5. Summary and conclusions. 

The FP of cluster galaxies is one of the most promising 
tools to study galaxy evolution, to constrain the epoch of 
galaxy formation, and to set constraints on the geometry 
of the universe. To those aims the FP has been recently 
studied in intermediat e-reds hift {z — 0.18 — 0.8 3) clusters 
(van Dokkum & Franx |1996| , Kelso n ct al. |l997| . Ben der et 
al. |1994 van Dokkum et al. |l998| , J0rgensen et al. |1999 , 
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Kelson et al. 200C). The main step in the above appli- 
cations is the comparison between FP determinations in 
different clusters, wavebands and redshifts, for which the 
determination of the FP coefficients and of the relative 
uncertainties plays a crucial role. 

The present study aims at clarifying various aspects in 
the determination of the FP coefficients and uncertainties. 



are systematic smaller than the true values and, what is 
more troublesome, are affected by a large scatter. 

The bootstrap technique gives, on the average, more 
accurate estimates than theoretical methods but for small 
samples the scatter of the individual values is much larger. 



Again, a minimum sample size of = Nm 
for the bootstrap to give reliable estimates. 



30 exists 



5.1. The problem of fitting the FP. 



5.3. Implications for the use of the FP. 



We introduced a statistical model of the FP that takes 
into account the measurement errors on the variables and 
the intrinsic scatter. The correlation between measure- 
ment uncertainties is also accounted for. 

We derived the MIST^ {Measurement errors and In- 
trinsic Scatter Three dimensional) fitting procedures, 
where i = logro, /Ze, logcio is the dependent variable of 
the fit. A bisector method, MISTbis , is also introduced. 

The assumptions under which the MIST procedures 
give unbiased estimates of the FP coefficients are that the 
errors and the intrinsic scatter do not depend on the 'lo- 
cation' on the fitting plane (HI), their covariance matrices 
are known (H2 and H3), and the covariance matrices of 
the observed quantities arc known (H4). 

Under these assumptions, the fits give the same coef- 
ficients whatever the choice of the dependent variable. 

For what concerns the above assumptions, the knowl- 
edge of the average measurement errors can be easily in- 
troduced in the procedures in order to satisfy HI and H2, 
while the problems posed by H4 can be overcome by as- 
suming a normal distribution for the observed quantities. 

Assumption H3, on the other hand, cannot be satisfied 
since we do not have a model for the intrinsic dispersion 
of the FP variables. We showed that it is precisely the 
lack of such a model that biases the best-fit coefficients in 
such a way that their values depend on the choice of the 
dependent variable and could not coincide with the slopes 
of the 'true relation'. 

5.2. Uncertainties on the FP coefficients. 

We addressed the problem of the estimate of the uncer- 
tainties on the FP coefficients by numerical simulations 
based on the best sampled FP so far (the Gunn r FP of 
Coma by JFK96). 

The results of the simulations are presented in Fig. 4 
in terms of 5a/a and 5h/h versus log(A^). This figure can 
be used to state the sample size needed to achieve a given 
accuracy for the FP coefficients. 

On the basis of the simulations we tested the perfor- 
mances of the methods used to estimate the uncertainties: 
the theoretical formulae and the bootstrap technique. 

We showed that for N > Nmin ~ 30 the theoretical 
formulae can recover the variances of the FP coefhcients 
within about 5%, while for N < Nmin the uncertainties 



As an application of the present study we addressed in 
Sect. 4 the question of the 'universality' of the FP. The 
main obstacles are the limited size of the available (homo- 
geneous) samples and their different completeness charac- 
teristics. 

We find that the slopes of the FP in the ten clusters 
considered are consistent within the large uncertainties 
(due to the small number of galaxies in most of the sam- 
ples). We stress that this result should be taken at most 
as an indication in favor of the universality of the FP and 
not as an evidence. 

With the available data samples it is not possible to 
settle the question of the universality: to this aim, larger 
and more homogeneous samples would be needed. 

For the application of the FP to the study of galaxy 
evolution it would be crucial to know the behavior of the 
slopes as a function of redshift. Pahre et al. (1998) built 
a model for the FP that predicts a change in slope with 
redshift: if age contributes to the tilt of the FP, massive 
galaxies should be much older than the least luminous 
ones. As a consequence of this 'differential' evolution, the 
slopes will change with redshift. They predict the change 
with redshift of the coefficient a of the FP (see their Fig. 
10) by adopting two different models for the stellar evolu- 
tion (Bruzual & Chariot 1998 and Vazdekis et al. 1996). 



The Bruzual & Chariot models produces a rapid change of 
a, that should decrease by about 0.15 from z=0 to z^O.3, 
namely, 12% of its local value in the Gunn r band. 

Recently, Kelson et al. ( ^000 ) derived the FP for a 
cluster at z=0.33 with 30 galaxies. The authors claim that 
the slopes are fully consistent with those of the local FP 
(Coma). They derive the FP by the orthogonal fit and 
give the following uncertainties for the slopes: 6a/a=0.1 
and i56/6=0.11, that are consistent with the values ob- 
tained here from the simulations with n=30 (see Fig. 4). 
These data seem thus to rule out the rapid evolution of a 
predicted with the Bruzual & Chariot model. 

The situation is quite different with the Vazdekis et al. 
models. At z — 0.3 the slope a decreases by only 0.03%. 
The minimum sample size needed to achieve such a pre- 
cision at z = 0.3 would be A^ ~ 300, which is out of the 
reach of any possible observation. The situation improves 
by moving at higher redshifts, where the magnitude of the 
predicted change of the slopes is larger. For the cluster at 
z=0.83 studied by van Dokkum et al. (|1998|) , in fact, the 
situation is by far more promising. To test the predicted 
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change in slope it would be sufficient to determine the FP 
on iV 40 galaxies, that is well within the possibility of 
a S-lOm class telescope in a reasonable number of nights. 

Appendix A 

From Eq. we write the probability density (hereafter 
PD) of Yi as: 



The symbols A, Ai and A2 will indicate the quantities 
obtained by substituting in Eqs. (B3 - B5) with un- 
biased estimates Cy. 

From Eqs. (O - Oh we have: 



Al _ Al 
A A 



fM)dY, = / fxM ■ ~ X,;X,)clX, 



(Al) 

where fviiVi) and fxi{xi) are the PDs of Yi and Xi, and 
[ipi] Xi)df>i is the probability that, for a fixed output of 
Xi, ^i is included between (jji and (f>i + dpi. For brevity, 

we put + = 

By Eq. (Al) we obtain the following relations: 



- Al -a A- A 



0}- 



13-13 ^ ^{(a2-A2) -/3(a-a)} + o 



(B6) 
(B7) 



= J fx^X,) j Y,U^{Y,~X,-X,)(¥, 



where a and (3 are defined by replacing Gij with Cy in 
Eqs. (0 , [n]), and o is a quantity that can be neglected 
when the sample size increases. 
By using the relation 



XrfxAXr) 



dX,, = 



E(I]) =E(S) 



(B8) 



fxAXr) 



f^i{4>i]Xi 



4'if<i>i{4>i;Xi)dj)i 



dX,, 



dX, 



(A2) 



V(r,) = / {Y^-l^{Y,f)fYAV^)dy^ 



where S is a generic RV and E(E) an unbiased estimate 
of E(E), the following relations can be written: 



Cy - Cy = E(r,r,) - E(y,y,) - E(y,)E(r,) + 

- E(r,)E(yj) + o (B9) 



ij^ki — '--'ij^ki — '-^ij I ^(YkYi) — Fi{YkYi) 



Cii c 



x^-E{X,y fxM 



/$ . {(j)i ; Xi 



dx, + + E(rfe) \^(Yi) - E(yo] - E(yo \^(Yk) - MYu)] } 



+ 2 / XJxAXr) 



+ / fxM 



(piUA'Pi'iXi 



'S'.Wi Xi)d(j)t 



dX, 



dX, 



(A3) 



Cki { [MY^Y,) - E(y,y,)] - E(y,) 
E(y,) [E(y,)-E(y,)]} + o 



E(y,)-E(y,) 



(BIO) 



By using the normalization of and the hypothe- 
sis E($,) = 0, from (A2) we obtain Eq. @). Eq. (|) fol- 
lows from (A3) with the further hypothesis that $i and 
{Xi} are not correlated ('hypothesis Tii'). 

If Til does not hold, we obtain again Eq. (||) pro- 
vided that the CM terms of {^mi} and {$si} are substi- 
tuted with their expected values with respect to the PDs 
of {X,}. 

In a similar way, Eq. (^ can be proved. 
Appendix B 

Let's start by introducing the following notations: 

= C($„,.,$„,J+C($s.,$s,) (Bl) 

= c(y„y,) (B2) 

A = [C12 - V12]' - [Cii - Vii] [C22 - V22] (B3) 
Al = [C23 — V23] [C12 — V12] — [Ci3 — Via] [C22 — V22] 

(B4) 

A2 = [C12 - V12] [Ci3 - V13] - [C23 - V23] [Cii - Vii] 

(B5) 



Setting Pijik = {Cij -Vij)YiYk and q^/ = 
(Cij — Vij)E{Yi), we introduce the RVs (, Ci and 
C2 by the following definitions: 

C ■ ^ = 2P1212 — P1122 — P22II + 2 [q221 ~ qi22] Yi + 

+ 2 [qii2 - qi2i] Y2 (Bll) 
Ci ■ ^ = P2312 + P1223 ~ P2213 ^ P1322 + [2qi32 — qi23+ 

- q23l] Y2 -I- [q223 - q232] Yi + [q221 - qi22] Y3 (B12) 

C2 • A = P1312 + P1213 ^ P1123 ^ P2311 + [2q231 ~ qi23+ 

- qi32] Yi -f [qii3 - qi3i] Y2 + [qii2 - qi2i] Y3 (B13) 

Using Eqs. (B8 - B13) and setting Cq = Ci ^ C/3 = C2 — 
/3C, after lengthy calculations, Eqs. (B6,B7) are rewritten: 



a~a ^ E(Ca) - E(Ca) + o 
P-13 ^ E(C/3)-E(C/3) + o 



(B14) 
(B15) 



From the above relations and Eq. (|T^), we also have: 
7-7 - E(C^)-E(C^) + o (B16) 
where = Y3 - aYi - PY2 - E(yi)Ca - E(y2)C/3. 



16 



La Barbera et al.: Fitting the Fundamental Plane 



If {Yi} are normally distributed, estimates of the ex- 
pected values in the previous equations are given by sam- 
ple means. The CM terms of a, (3 and 7 are thus obtained 
by the usual formulae: 

" NiN- n ^ ^^^^ ~ - 6) 

' i=l,N 

^ a,f3,'j (B17) 

where N is the sample size and Cq, Cp and are the RVs 
obtained by substituting in the expressions of (^c, (^^ and 
C7 the unknown quantities with their unbiased estimates. 
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